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Stellar dynamical evidence against a cold disc origin for stars in the 
Galactic Centre 
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ABSTRACT 

Observations of massive stars within the central parsec of the Galaxy show that, while most 
stars orbit within a well-defined disc, a significant fraction have large eccentricities and / or 
inclinations with respect to the disc plane. Here, we investigate whether this dynamically hot 
component could have arisen via scattering from an initially cold disc - the expected ini- 
tial condition if the stars formed from the fragmentation of an accretion disc. Using N-body 
methods, we evolve a variety of flat, cold, stellar systems, and study the effects of initial disc 
eccentricity, primordial binaries, very massive stars and intermediate mass black holes. We 
find, consistent with previous results, that a circular disc does not become eccentric enough 
unless there is a significant population of undetected 100-1000 M objects. However, since 
fragmentation of an eccentric disc can readily yield eccentric stellar orbits, the strongest con- 
straints come from inclinations. We show that none of our initial conditions yield the observed 
large inclinations, regardless of the initial disc eccentricity or the presence of massive objects. 
These results imply that the orbits of the young massive stars in the Galactic Centre are largely 
primordial, and that the stars are unlikely to have formed as a dynamically cold disc. 
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1 INTRODUCTION 

Observations of the Galactic Centre (GC) have revealed a popula- 
tion of young (~5 Myr) massive stars orbiting within a fraction of 
a parsec of the M BH ~ 3 x 10 6 M G supermassive black hole Sgr A*. 
The existence of any young stars so close to the black hole is a 
puzzle, since the tidal field would prevent star formation from gas 
clouds with densities p < M Bli /[(4n/3)r 3 ]. This critical density is 
~ 10 10 cirT 3 at 0. 1 pc, far higher than the densities of normal molec- 
ular clouds. The dynamics of the stars are also unusual and in- 
teresting. Although the innermost 'S-stars' appear to be consistent 
with isotropic orbits, many of the yo ung stars at a distance of the 
order of 0. 1 pc lie on a single plane JLevin & Beloborodov||2003j; 
iGenzel et alll2003l : IPaumard et al.ll2006l; |Lu" et al. 2006]) and rotate 



as e ~ 0.8 for stars on counter-clockwise orbits ( IPaumard et al.l 
2006Q 

Two classes of theories have been advanced to explain the ori- 
gin of the disc stars. In the first, the stars form in situ from a mas- 
sive accretion disc in which the d ensity is large enou gh to exceed 
the tidal limit JLevin & Beloborodov 2003; Milosavljevic & Loeb 
l2004l ; lNavakshin & Cu adra 2005; Paumard et al. 2006]). Accretion 
discs at these radii around black holes are known to be vulnerable to 
fragmentation and star formation provided only that they can cool 



on a dynamical timescale (Paczyriski 1978; Kolykhalov & Sunyaev 
1980; Shlosman & Begelmanl ll989J : iGammiel l200ll : iRice et al 



clockwise on the sky. The eccentricities of stars in this clockwise 
disc are significant. The average value is e » 0.35, with the eccen- 
tricities of sev eral stars going above 0.5 1 B eloborodov et al.1120061 : 
iLu et alj|2006l) . Moreover, the stars that do not belong to the clock- 
wise disc (which, by definition, have large inclinations that in many 
cases exceed 90°) have even larger eccentricities, perhaps as large 



120031 ; iGoodmanl 120031) . The second model posits that the stars 
formed outs ide the GC region but migrated in due to dynamical 
friction (e.g.lGerhardlfe oOl; Ki m & Morrisl|2003l; JKirri et alj|2004 
iGtirkan & R asio 2005; Berukoff & Hansen 2006). Understanding 
the predicted dynamics of these models in more detail may assist 
in discriminating between them. 

Fragmentation of a thin accretion disc yields, almost by con- 
struction, starts with the basic disc-like dynamics observed in the 
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1 Note that there is some dispute as to whether there is evidence for a sec- 
ond, counter-clockwise disc of stars. This question is not important for our 
purposes. We rely on observations only insofar as they demonstrate that 
there is (at least) one disc, and that some stars have orbits that do not coin- 
cide with any suggested disc plane. These aspects are uncontested. 
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GC. Within such a model those stellar orbits that do not belong 
to any well-defined disc have to be generated after the disc frag- 
ments, perhaps via scattering. There is, however, only just enough 
time f or significant dynamical evolution to occur. I Alexander et al.l 
(2007) found that even with a top-heavy mass function it was 
hard to attain the large eccentricities observed. This problem can 
be circumvented if the disc fragments while it is still eccentric 
JNavakshin et~aT1l2007l : Alexander et"alll2008h . but it is still diffi- 
cult to explain the high inclinations. Here, we investigate whether 
stellar dynamical evolution from any plausible cold disc initial con- 
ditions (i.e. those for which the local velocity dispersion is small 
compared to the orbital speed) can explain the observed stars. We 
consider both circular and eccentric initial conditions, and allow 
for the possibility of a high primordial binary population and addi- 
tional scattering from unobserved massive objects (short-lived mas- 
sive stars, or intermediate mass black holes). 



2 NUMERICAL METHOD 

We model the stellar dynamics of a disc of stars using the mercury6 
code, which is optimiz ed to study systems dominated by a central 
mass dchambersl 1 1 9990 . In order to treat close encounters reliably 
we adopt the optional Bulirsch-Stoer integrator, in its version for 
conservative systems. 

We simulate stellar discs extending radially from 0.05 to 
0.15 pc. |3 The discs are composed of 100 stars randomly dis- 



sinqle stars 



tributed following a surface density profile E(r) 



roughly 



as observed in this radial range, with the mass of each star equal 
to 25 M = 8.3 x 10~ 6 M BH . In most of our models the stars are 
initially rotating in circular Keplerian co-planar orbits around the 
central object, with small random perturbations such that the thick- 
ness of the disc is h/r = 0.01. The exact value of this perturbation 
is not important, as the evolution of the disc is initially rapid. We 
integrate the systems for a time roughly equivalent to 2700 orbits at 
r = 0. 1 pc, almost 5 Myr. Since N is rather small, we run three sim- 
ulations with each set-up to monitor the variance in the outcome. 

Our model system includes the number and masses of the ob- 
served population of stars in the GC. Alth ough observational ev- 
idence favours a top-heavy mass function tNavakshin & Sunvaevl 
120051 ; IPaumard et~al] |2006) there could still be undetected lower- 
mass stars in the same region. Low mass stars, however, al- 
ways act to damp the dispersion of the observed massive stars 
JAlexander et alJ uOOT). so in ignoring them our models yield an 
upper limit to the rms eccentricities and inclinations that result from 
scattering. 



3 SIMULATIONS 

3.1 Single stars 

To establish a baseline for subsequent comparisons, we first model 
the evoluti on of a disc of single stars similar to that considered 
previously ( [Alexander et al.l 120071) . We also use this simple sys- 
tem to test the numerical convergence of the integrator, by run- 
ning three sets of models with different integrator tolerances. Each 
set includes three independent realisations of the initial conditions. 



2 Young stars are observed further away, but their density decreases rapidly 
and their scattering time-scales become longer with radius, so we do not 
expect the omission of stars at r > 0. 15 pc to alter our results. 




Figure 1. Eccentricity (upper curves) and inclination (lower curves) rms 
values from the nine realisations of discs of single stars. There is a small 
bias toward larger dispersion in the runs with less stringent accuracy re- 
quirements. For comparison, a curve e oc r 1 ' 4 is also shown. 



Figure Q] shows the evolution of the rms values of the eccentric- 
ity and inclination for each simulation. We find that convergence 
is obtained for an accuracy parameter specified by accuracy = 
2 x 10~ 13 , and that results accurate to the 10% level, adequate for 
our purposes, can be obtained with accuracy = 10~ 12 . We use the 
latter value for most of our subsequent runs. 

As expected theoretically, the eccentricity and inclination rms 
grow roughly as t l/4 , with a proportionality e rms » 2/ rms (e.g., 
iLissauenl 19930 . where i is measured in radians. To compare quan- 
titatively with previous work, we fit the average eccentricity rms 
from all nine simulations with the analytical solution, 



dcr, 
IT 
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where a, is the velocity dispersion of a stellar disc of N, stars 
each of mass M„, orbiting at radius r with radial extension Ar. The 
quantity In A, is the Coulomb logarithm, t olb the o rbital period at r, 
and C i a geometrical constant of order unity (e.g jAlexander et al.l 
120070 . We find a value of Ci = 1.65, in good agreement with the 
1.85 that I Alexander et alj ( 120070 found from their high resolution 
runs. This is a very small difference in practice, as the eccentricity 
rms value scales only as C7 . In agreement with prior work, we 
conclude that an initially circular disc set up to match the numbers 
and masses of observed stars in the GC attains an rms eccentric- 
ity of around 0.2 within a few million years. This is too small to 
explain the larger eccentrities observed in the GC population. 



3.2 Primordial binaries 

A small nu mber of bin ary stars have been identified in the Galac- 
tic Ce ntre JMartins et a l. 2006; Peepl es et aT]|200 7; Rafels ki et al.l 
120070 , and, given our current knowledge of disc fragmentation, it 
is possible that the primordial binary fraction resulting from a frag- 
menting disc could be significant. Encounters between stars and 
hard binaries tend to shrink the binary orbit, releasing energy that 
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Figure 2. Eccentricity rms values from the nine different simulations of 
discs of binary stars. The tighter the binaries are initially, the more disper- 
sion the system acquires. 



Figure 3. Eccentricity (upper curves) and inclination (lower curves) rms 
values from the eccentric disc runs. The eccentricity rms is roughly constant 
and the inclinations remain small in the more realistic case (see Section |33V 



can heat the system dHeggie|[l975l : IPerets euill feoOSV 3 ! To study 
whether binary heating could be significant in the GC we consider 
an idealised model in which the initial binary fraction is 100%. We 
replace the disc of 100 single stars with one made up of 50 equal 
mass binaries, whose centres of mass are initially in circular orbits 
about the black hole. For the fiducial case the binary orientations 
lie in the disc plane. In each simulation all binaries have the same 
initial value of semi-major axis, which we vary between 10~ 3 pc 
and 10~ 4 pc. A binary with semi-major axis 10~ 4 pc is hard (and 
hence will shrink and give energy back to the system) if the stellar 
velocity dispersion rx, i 20kms _1 , which corresponds to e ~ 0.1. 
For numerical reasons we are only able to track binaries as they 
shrink down to a minimum size which we set at 3 x 10~ 6 pc. If any 
stars get closer than this limit, the code merges them into a single 
object conserving energy and momentum. 

Figure|2]shows the evolution of the rms value of the eccentric- 
ity for the binary runsQ When the binaries' initial separations are 
10~ 3 pc they are only marginally bound and most of them quickly 
get disrupted. As a result the system behaves as if it were made up 
of single stars from the start, and the resulting eccentricities and 
inclinations are essentially identical to those discussed in the pro- 
ceeding section. For initially tighter binaries, we see a trend where 
closer binaries produce more eccentric orbits. The plots also get 
noisier, as the instantaneous Keplerian orbital elements of a few 
stars are strongly perturbed during close encounters. The changes 
are relatively modest, and even a disc made up entirely of the tight- 
est binaries yields a final eccentricity that is only ~ 25% larger 
than the single star case. The rms value of the inclination is again 
roughly half of the eccentricity. Initial binary orientations that are 
not aligned with the disc plane yield modestly smaller inclinations 
and eccentricities than the aligned case. We interpret this as im- 
plying a reduced cross-section for binary-binary encounters in the 



3 Supernovae can also heat the disc, as t he s urviving companion keeps the 
velocity it had in the binarv. lPerets et alj 120081) , however, showed that even 
with 100 km s kicks no significant inclinations are reached. 

4 The average here is defined not over stars, but rather over 'objects' that 
could be either bound binaries or single stars. It is very hard to determine 
whether the observed stars are actually singles or binaries, so the value we 
are plotting is what can be compared to current observations. 



non-aligned case. We also experimented with decreasing the min- 
imum distance below which stars are assumed to merge. A reduc- 
tion in this scale by a factor of three resulted in slightly smaller 
eccentricities, but we were unable to demonstrate convergence by 
evolving systems with even smaller merger radii due to numerical 
difficulties. Accordingly, we do not attach physical significance to 
this difference. In reality, the dynamics of any binaries that shrink 
to such small scales (which are less than an AU) would likely be 
affected by the presence of circumstellar discs left over from their 
formation. We cannot model such effects, but it is hard to see how 
they could increase the overall stellar disc eccentricity. 

3.3 Eccentric flat discs 

If the stars originated from the fragmentation of a self-gravitating 
disc, they did not necessarily form after the disc became circular. It 
is possible that gas settled in a plane around the black hole and frag- 
mented before it had time to circularise. In this case, the stars would 
initia l ly have inherited the e ccentricity of the gas ( INavakshin et alj 
120011 Alexander et alj|2008h . 

To study the growth of inclinations in this scenario, we ran 
several simulations with initially eccentric stellar discs. The stars 
were set up in orbits with e = 0.5, with a surface density E(a) oc a 
between 0.05 and 0.15 pc. We considered both the case where 
the angle of pericentre of the eccentric orbits was initially con- 
stant (so all the ellipses are nested), and the case where this angle 
was randomly distributed (in this case the in-plane velocity disper- 
sion is initially large). The aligned case is the expected outcome 
immediately after the fragmentation of an eccentric gaseous disc, 
but the M CUSD ~ 10 5 M Q worth of old stars between 0.05-0.15 pc 
JGenzel et alj 12003!) will lead to differential precession and ran- 
domly oriented orbits on a time-scale ~ (Mbh/^cusp)^*- This time- 
scale is very short compared to the current age of the stellar system, 
so the more realistic set up for our simulations without the cusp po- 
tential are those with randomly chosen pericentre angles. 

Figure [3] shows the results for these simulations, with the blue 
curves depicting the case of randomly aligned orbits. The typical 
eccentricity remains fixed at e « 0.5, while the inclination reaches 
a typical value of 0.1, just as in the circular case (Sect. [3~TV The 
runs starting with aligned orbits (red curves) have a larger spread 
in eccentricities (although still e rms as 0.5), with approximately 
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Figure 4. Eccentricity (upper curves) and inclination (lower curves) rms 
values from the six different simulations of discs of stars that included a 
significant number of very massive stars. The very massive stars, however, 
are not included in this analysis, as they are not observed and are supposed 
to have already died. The typical eccentricity reaches an interesting value, 
but the inclinations are still too low to match the GC observations. 



Figure 5. Aitoff projection of the angular momenta of the 25- M Q stars in 
a simulation that also included 100 100- M Q stars. The angle definition is 
chosen such that the disc initially has the same orientation as the observed 
clockwise disc (cf., Lu et al. 2006, Fig. 5). No star reaches the bottom half 
of the plot, where they would be observed rotating counter-clockwise. The 
inset shows the distribution of semi-major axes and eccentricities. 



10 of the innermost stars acquiring large eccentricities and incli- 
nations as a result of the interaction with the eccentric potential 
of th e rest of the stars. This se ems to be a resonant relaxation ef- 
fect (R auch & Trema indll996T) , a mechanism already invoked by 
iHopman & Alexa nder (2006) and lLevinl J2007i) to explain the high 
eccentricities and inclinations of the innermost 'S-stars'. The frac- 
tion of the stars that get high i orbits is still too small to match the 
GC population, and, as explained above, we expect the background 
stellar potential to substantially weaken this resonant process. 



3.4 Influence of Massive Stars 

So far we have concentrated on the stars that are currently observed. 
However, it is possible that many more massive stars were born 
but have already exploded as supemovae. To take into account this 
effect, we ran a set of simulations where - in addition to the 100 
25- M stars previously considered - we include an extra number 
of 100 M G stars. These stars are included throughout the simulation 
but are not taken into account in the final analysis, so in practice it 
is as if these very massive stars disap peared a very short time ago. 
As shown by I Alexander et al.l ( 120071) . the more massive stars will 
sink toward the mid-plane, while the lower mass stars will acquire 
higher eccentricities and inclinations. 

We ran several simulations with either 50 or 100 of the 
100 M Q stars. Figure [4] shows the eccentricity and inclination rms 
for the 25 M Q stars. The eccentricities now reach a value e w 0.4, 
comparable to that observed. The inclinations, however, are still 
too low to explain the non-disc population in the Galactic centre. 
There is not much difference between the cases with 50 and 100 
very massive stars, and, although we cannot provide a strict limit, it 
is hard to imagine that even more massive stars were formed in the 
GC. It does not seem likely that a significant non-disc population 
can be formed via scattering off such objects. 

To illustrate the failure of the model more clearly, Figure [5] 
shows the orbital parameters and angular momentum direction of 
each 25 M G star from a run with 100 extra 100 M Q stars. Even in 
this case only a handful of stars seem to be outside of the disc, and 



none of them appears to be rotating in the opposite direction. The 
eccentricities, however, are comparable to the observed values. 



3.5 Intermediate Mass Black Holes 

Finally, we investigated the possibility that scattering from a hypo- 
thetical po pulation of intermedia t e mas s black holes might heat the 
disc stars. IPortegies Zwart et alj ( 120061) argued that core collapse 
in GC star clusters could yield a significant population of such ob- 
jects, with masses around 1000 M Q . Adopting their numbers, there 
could be ~50 intermediate black holes within 10 pc of the GC, and 
some would likely have orbits that intrude upon the observed stellar 
disc. 

To study the effect of these black holes, we set up simula- 
tions in which the stellar disc was initialised as described in Sec- 
tion l3.ll We then added N intermediate mass black holes, with or- 
bits appropriate for an isotropic cluste r with number density pro- 
file n(r) oc r~ 1/4 ( IBahcall & Wolf 1976) (i.e. drawn from an energy 
distribution f(E) oc E l/4 with random directions). We consider a 
model with N = 5 1000-M a inside 1 pc (similar to the number 
expected from lPortegies Zwart et alj|200q analysis), together with 
additional cases that assume N = 10 500-M o and N = 20 500-M o 
intermediate mass black holes. 

The resulting £„ curves are shown in Figure [6] There is a 
lot of dispersion between the different realisations, since the small 
numbers of black holes involved mean that in some cases none (or 
very few) of the massive objects intersect the stellar disc. However, 
it is clear that the typical effect is rather similar to the case of mas- 
sive stars discussed previously. Quite large eccentricities are readily 
produced, but the inclinations (i lms » e !ms /2, not shown in the plot 
for clarity) remain too small to explain the observed distribution of 
stars in the Galactic Centrqj. 



5 Note that although iKozail j 19621) resonance would allow for the gener- 
ation of high inclinations, the rapid precession of stellar orbits due to the 
stellar cusp is likely to destroy this resonance. 
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Figure 6. Eccentricity rms values from the nine different simulations of 
discs of single stars perturbed by IMBHs of different masses. The results are 
similar to those of the models with extra massive stars in the disc (Fig. [4). 



Ideas that would yield dyna mically hot initial conditions - such as 
very strongly warped discs (Nayakshin 2005) or clouds that col - 
lapse as they engulf the black hole ( Yusef-Zadeh & Wardle 2008), 
perhaps forming two coeval discs (Genzel et al. 2003) - have been 
proposed, though whether such scenarios can be analysed in the 
same framework as the simpler disc models is dubious. 

It is to be hoped that better observational data on the stellar 
orbits, which should confirm or refute the existence of additional 
clustering in the distribution of orbital vectors, will afford some 
insight into which models are viable. 
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4 CONCLUSIONS 

In this paper we have studied the dynamical evolution of initially 
flat cold stellar discs consisting of the massive stars observed today 
in the GC. We find that these systems do not evolve to match the 
observed stellar population - which includes significant numbers of 
stars that are not clearly identified with a well-defined disc plane - 
even when scattering from plausible (and perhaps not so plausible) 
additional perturbers is included. It is clear that the most important 
observational constraints on theoretical models comes from mea- 
surements of the orbital inclinations. Although a flat circular disc 
does not become eccentric enough in 5 Myr to match the observa- 
tions, it is possible to create quite eccentric systems either by start- 
ing from eccentric gas discs, or by perturbing circular discs with 
additional massive objects (short-lived stars, or intermediate mass 
black holes). However, in no cases that we have considered is it pos- 
sible to generate enough very high inclination or counter-rotating 
stars within the available time. No physical process we have con- 
sidered breaks the basic relation e lms a 2z lms between the mean 
eccentricity and inclination and this, coupled with the slow rate of 
heating at late times (typically oc ? 1/4 ), precludes formation of a 
significant non-disc population. 

The simplest disc scenario for forming the sub-pc stars in the 
GC appeals to the fragmentation of a flat gaseous disc into a pop- 
ulation of massive stars. Both the tendency of self-gravitating ac- 
cretion discs at these radii to fragment, and the subsequent forma- 
tion of a top- heavy mass function, receive a t least qualified support 
from theory (Shlosman & Begel man 1 19891 : iLevin & Beloborodovl 
12003c iNavakshin et alj|2007l) . However, our results suggest that it 
is very unlikely that subsequent evolution of a stellar system formed 
in this manner would yield any substantial population of stars 
outside of the original disc plane or planes. Possible alternatives 
are that some or all of the GC stars migrated inward f rom larger 
radii (though this model has its o wn challenges, e.g., iKim et al.l 
12004c INavakshin & Sunvaevl 120051) : that a non-spherical compo- 
nent in the background stellar population, or a second stellar disc 
( INavakshin et al.l 2006), influences the dynamics of the main stel- 
lar disc in ways we have not considered; or that in situ star for- 
mation occurs in a very complex geometry. In the latter case the 
present orbits of stars are expected to trace, to a first approxima- 
tion, the angular momentum distribution of the star forming gas. 
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